Limits on the amplification of evanescent waves of left-handed materials 
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We investigate the transfer function of the discretized perfect lens in finite-difference time-domain (FDTD) 
and transfer matrix (TMM) simulations; the latter allow to eliminate the problems associated with the explicit 
time dependence in FDTD simulations. We argue that the peak observed in the FDTD transfer function near 
the maximum parallel momentum fc|| jmax is due to finite time artifacts. We also find the finite discretization 
mesh acts like imaginary deviations from fj, = e = — 1 and leads to a cross-over in the transfer function from 
Constance to exponential decay around ft|| jmax limiting the attainable super-resolution. We propose a simple 
qualitative model to describe the impact of the discretization. ftii max is found to depend logarithmically on the 
mesh constant in qualitative agreement with the TMM simulations. 
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The ability of the left-handed finite slab with a homoge- 
neous permeability fi = — 1 and permittivity e = — 1 to form 
a perfect lens (PL) has received much attention since first pro- 
posed by Pendry 1 . Such a slab does not only compensate 
the phase of the propagating waves emanating from a point 
source to form a focus on the opposite side of the slab. It also 
amplifies the evanescent waves which decay exponentially in 
vacuum into exponentially growing solutions inside the slab. 
This way all the source amplitudes reemerge in the focus. The 
immediate consequence of this behavior is that the resolution 
of the image may overcome the diffraction limit. Soon after, 
it was realized 2 ' 3 that the restoration of the evanescent waves 
by the PL is exceptionally sensitive to small deviations from 
fi = e = — 1. The transfer function of the PL, defined by 
the amplitude ratio of a plane wave component at the focus 
and the source, is, in the ideal case, unity for all to and ku 
up to infinity. For the near-perfect lens it exposes an order 
of unity (o(l)) behavior at small parallel momenta ku which 
turns into exponential decay ~ e~ k u d for large ku. The cross- 
over between o(l) behavior and exponential decay for a given 
PL defines a maximum parallel momentum fcii max which 
qualitatively constitutes the highest evanescent wave still re- 
stored by the PL, hence, defines the maximum attainable sub- 
wavelength resolution Ax m i n ~ 2ir/kn max.- For small devia- 
tions fj, = e = —1 + 7 with 7 £ C, \j\ <C 1 from the ideal PL 
a logarithmic dependence kn max d ~ — log |7| of the cross- 
over momentum has been found 3 ". Here and throughout the 
paper we employ a dimensionless formulation measuring all 
lengths in units of the linear size L of the unit cell and all fre- 
quencies in units of the vacuum speed of light divided by L. 
In particular, this renders the dimensionless vacuum speed of 
light c = 1 and wavelength A = 2ir/u. 

Almost all numerical investigations of the PL's imag- 
ing properties deploy finite-difference time domain (FDTD) 
simulation using a time and space discretized version 
of the Maxwell equations. After a few contradictory 
publications 4 - 5 - 6 , Rao and Ong 7 - 8 established the amplifica- 
tion of the evanescent waves inside the LH material slab and 
the cross-over behavior in the transfer function numerically 
for the FDTD method. They also observed the occurrence 



of surface plasmons, ie. local field enhancement, at the first 
interface for the slightly lossy PL. The FDTD simulations 
of the PL suffer from explicit time dependence. Since the 
time domain simulations involve a finite time window from 
the "switch-on" to the actual measurement of the fields, the 
results are obtained as a superposition of a finite width u>- 
distribution around the target frequency luq which narrows 
as the simulation time increases. The corresponding transfer 
function £fdtd (^o) differs considerably from the stationary 
(frequency-domain for a single frequency u>) transfer function 
too(w). In conjunction with the physically always present dis- 
persion fJ,(u>), e(ui) of the left-handed material this leads to a 
possible coupling to the surface plasmons on both interfaces 
of the slab which in turn causes convergence problems in the 
FDTD. The FDTD only converges for the near-perfect lens 
where the surface plasmons are damped by the small imagi- 
nary part in the LH slab 8 -. If we approach the ideal PL the 
FDTD ceases to converge which renders the method unusable. 

Due to the existence of surface plasmons at the LH slab 
the transfer function of the near-perfect lens includes poles 
along the surface plasmon dispersion relation 9 which can be 
approximated by — tanb.(fei|<i/2) = — l±7forfeii well above 
the propagating modes. By virtue of the LH materials disper- 
sion relation this essentially real 7 directly translates into a 
frequency deviation from ujq. The poles approach loq expo- 
nentially for growing fcii . For small ku the poles of the surface 
plasmons are usually outside the finite width tj-distribution, 
we find convergence of the FDTD and the transfer function 
is dominated by the stationary transfer function too(u)) at 
ujq. For large k\\ the poles are exponentially damped in all 
lossy cases and cease to contribute either. However, for in- 
termediate ku « fcii max the surface plasmon poles constitute 
the principal contribution to £fdtd(wo)- This leads to non- 
convergence of the FDTD due to the emerging "beating pat- 
tern", modulated by the frequency difference of the two sur- 
face plasmon branches as explained by Gomez-Santos^. As 
we emphasize here, this also explains the unexpected peak 
iFDTD(wo) 3> 1 around fc|| max in the FDTD transfer function 
found by Rao and Ong 7 and also confirmed by our own FDTD 
simulations. The peak originates from the contribution of the 
diverging t 00 (w po i c ) to the magnitude of iFDTD(wo). Ana- 
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lytically we would expect a monotonous transition from the 
o(l) behavior below the cross-over momentum below fc|| max 
to exponential decay above for the near-perfect lens. 

In general we are interested in the stationary case transfer 
function too(k\\,ui,d) of the PL as this allows us to estimate 
the imaging properties. The field components in the focus are 
given by the field components in the source as 
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where g(uj, uiq) is the frequency distribution around the fre- 
quency of the point source due to the switch-on of the source 
and the finite observation time window. E SOUIce (k\\) and 
g(LU,LU ) are parameters of the setup and measurement, only 
toc(k\\,oj) is an intrinsic property of the lens. In the FDTD 
this stationary transfer function too(k\\,u}) is only accessible 
via the temporal Fourier transform of the simulation results 
which is neither convenient nor especially robust against nu- 
merical error. 

In contrast, the transfer matrix method (TMM) simulations 
provide a means to directly obtain the stationary transfer func- 
tion for a single frequency. We can simulate transmission and 
reflection amplitudes as well as the spatial field distribution 
without the artifacts of finite time simulations. In order to 
eliminate the effect of numerical error we employ an arbitrary- 
precision implementation of the TMM described for instance 
in Ref.[I3 

Fig. [Qa) shows the TMM simulated spatial field distribu- 
tion for the parallel component of the E field in TE polar- 
ization across a PL for one particular evanescent wave com- 
ponent with fen = 1.89 u> and several spatial discretizations 
ranging from 240 to 12857 linear mesh points per vacuum- 
wavelength. The LH slab extends from z = 3 to z = 9, which 
corresponds to a thickness of 0.28A. The vertical solid lines 
indicate the interfaces. The source is located at z — and we 
get an image at z = 12. Note that despite the elimination of 
the explicit time-dependence of the FDTD simulations we still 
observe unexpected field enhancement at the first interface of 
the PL. These artifacts are easily associated with the finite dis- 
cretization leading to reflection of incident evanescent modes 
at the first interface. Only for very fine discretization meshes 
the field distribution approaches the analytically expected zig- 
zag-form featuring a minimum at the first interface. For coarse 
discretizations we observe a prominent maximum at the first 
interface accompanied by (almost) zeros of the fields before 
and after the interface indicating a phase shift of the response 
of the interface. In Fig.^b) we show the corresponding field 
distributions for the lossy PL where a small imaginary part 
7 = 0.03 i is added to the permeability and permittivity of the 
LH slab. Again we observe the field enhancement at the first 
interface, but the zeros of the field have disappeared. How- 
ever, in this case the behavior is dominated by the losses in 
the LH slab and the dependence on the discretization is much 
weaker. The observed dependence on the imaginary part for 
fixed discretization (not shown) confirms previous results ob- 
tained from FDTD simulations 7 . 

For the lossy PL there is a simple physical explanation 
for the reflection of evanescent waves and the occurrence 
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FIG. 1: (Color online) The distribution of the electric field in TE 
mode is shown from source (z = 0) to focus (z = 12) across the 
lossless (a) and lossy (b) perfect lens for successively finer discretiza- 
tion. The left-handed slab between the interfaces at z — 3 and 2 = 9 
has fi — e — — 1 for the lossless PL and fj, — e — — 1 + 7 with 
7 = 0.03 i for the lossy PL. We discretized using a uniform cubic 
mesh with a linear resolution ranging from 240 to 12857 mesh points 
per vacuum-wavelength. 



of surface modes at the first interface: For the evanescent 
waves in vacuum between source and first interface fcii is 
real and k± purely imaginary. Whenever the PL involves 
an (causal) imaginary part or the evanescent solutions on 
the right hand side of the slab couple to propagating modes 
or are subject to absorption, electromagnetic field energy is 
dissipated in the system. This energy has to be provided 
by the source and transmitted across the vacuum gap before 
the PL. Although it is well known that a single evanescent 
wave cannot transmit energy this is not true for a superposi- 
tion ^4 e i ( k ii+ k i) r + 5 e *(k||-ki)r of i ncom i n g an d reflected 

evanescent wave component. The general equation for the 
time-averaged Poynting vector for the TE mode inside the 
vacuum slab is 
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and similarly for the TM mode. Considering the first term it 
is immediately clear that in order to have an energy current 
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FIG. 2: (Color online) The transfer function from source to focus is 
shown for the lossy PL (symbols) with /j, = e = — 1 + 7 for two dif- 
ferent small imaginary parts 7 = 0.002 i (a) and 7 = 0.005 i (b) and 
different discretizations. The dashed lines show the corresponding 
transfer function for the lossless PL. We discretized using a uniform 
cubic mesh with a linear resolution ranging from 643 to 51429 mesh 
points per vacuum-wavelength. The additional rightmost dashed line 
corresponds to 102858/ A. 



normal to the interface across the gap, lm(AB*) and thus the 
reflection at the first interface has to be non-zero. Note that 
this even applies to the ideal PL if we have an outgoing energy 
current or dissipation on the right hand side of the lens. 

Now we shall consider the transfer function t 00 (kn,ui,d) 
of the PL with and without losses from source to focus as 
obtained by the TMM. In Fig. |2 we show the dependence 
of the transfer function for a fixed frequency u = 3/10 on 
the parallel momentum feii for several spatial discretizations 
for the lossless PL and two lossy PLs with imaginary parts 
7 = 0.002 i and 7 = 0.005 i added to both the permeabil- 
ity and the permittivity of the LH slab. Let us first consider 
the lossless PL represented by the dashed lines in both panels. 
For all discretizations there is clear evidence for a cross-over 
from o(l) behavior to exponential decay in the kn dependence 
of the transfer function. The cross-over occurs monotonously 
without a peak near fcii roax and kn raax increases with finer 
spatial discretization. This indicates that the discretization 
mesh constant acts like an effective imaginary part in a contin- 
uous lossy PL. For the lossy discretized PL, ie. adding an ex- 
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FIG. 3: (Color online) The logarithmic scaling of the cross-over par- 
allel momentum fcii max with the number of linear mesh points per 
vacuum-wavelength is shown for the lossless (7 = 0) and two lossy 
PLs with 7 = 0.002 i and 7 = 0.005 i. The dotted line is a fit 

&||,maxd = (24/27) log (4A/JV mcsh ) for the lossless PL. 



plicit imaginary part 7 we observe the same qualitative behav- 
ior. Here the cross-over is determined by both, discretization 
and the losses due to the explicit imaginary parts. For small 7 
and coarse discretization the behavior of the transfer function 
is entirely dominated by the finite discretization: the lossy PL 
virtually coincides with the transfer function for the lossless 
PL. For successively finer discretizations 7 starts to dominate 
the behavior, leading to a saturation of the discretization de- 
pendence of kn max at a value determined by 7. These results 
show that for a given lossy PL there is always a minimum dis- 
cretization mesh constant where the transfer function "con- 
verges", ie. becomes independent on the discretization. For 
the simulated lossless perfect lens the cross-over in the trans- 
fer function due to the discretization becomes the primary lim- 
iting factor for the observation of sub-wavelength resolution. 
In Fig. [3] we show the dependence of the cross-over momen- 
tum ifeii max on the discretization for the lossless and two lossy 
0.28A-PLs at u> = 3/10 as extracted from the data presented 
in Fig. [2] It is evident that for the lossless PL over a wide 
range of discretizations fcii max increases logarithmically with 
the linear number of mesh points per vacuum wavelength. For 
the lossy cases this slow increase saturates at a finite feii )niax 
which in turn decreases with increasing deviation 7 from the 
lossless case. In order to achieve a moderate five-times bet- 
ter resolution than the one provided by the propagating modes 
alone for the d = 0.28A lossless PL, we have to push the 
discretization to a ridiculously high value of 10 5 linear mesh 
points per vacuum wavelength. Such discretization mesh den- 
sities are easily limited by the available computer power. 

The effect of the discretization can be qualitatively under- 
stood in terms of a simple model. In the standard discretiza- 
tion of the Maxwell equations the E and H field components 
are assigned to the links of two mutually dual lattices**. As a 
consequence a wave traveling towards the surface of a dis- 
cretized homogeneous slab will first "see" the electric re- 
sponse and approximately half a mesh step later the magnetic 
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response (or vice versa, depending on the material discretiza- 
tion and definition of the interface). This can be analytically 
modeled assuming a continuum lossless PL with /i = e = — 1 
to be sandwiched between two thin layers with \i = —e= 1 
and — \i = e = 1, respectively. The thickness 8 of the sur- 
face layers shall be of the order of the discretization mesh 
constant. Now we can derive the leading order <5-corrections 
to the transfer function analytically using the transfer matrix 
technique. We can calculate the total transfer matrix of the 
left-handed slab (t 2 ) wrapped in surface layers (ti,t 3 ) and 
the two surrounding vacuum slabs (r 2 ) as 



imaging 



r (b) [^(8)^)^(8)}^). 
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For the transfer matrix of an homogeneous slab in wave rep- 
resentation we find 



Ti(d) = 



ai {d) pi(-d) 
(3i(d) cti(-d) 



with the elements 

a t {d) 

Hi{d) 



COS(M) + o ( Ci + 



Ci 



sin(kid). 



sm(kid) (3) 
(4) 



The Ci are defined as Q = fJ,ik / {noki) or Q = e h/(eiko) 
for the TE and TM mode, respectively; indices zero refer to 
quantities in vacuum. The transfer function coincides with the 
transmission coefficient t (we choose t- for convenience) of 
the imaging scattering matrix S imaging = S [r imaging ] , 



[l + o(S)] e lko< - a 
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cos(fc 2 d) + l + <5 2 (^^) isin(k 2 d) 



(5) 



We immediately recognize that the surface correction 
8 2 uj a / (oj 2 — k 2 ) = 5 2 lo a /k 2 in the denominator acts like an 
imaginary part in the permeability or permittivity of the near 
perfect lens. If the perfect lens condition a + b = d is satisfied 



we have t_ 



(l + o(8))/(l + 8 2 co% 2 (l-e 



-2ik2 d 



)). Let 



us now consider the transfer function for evanescent waves, 
ie. fc|| > u. Then fc 2 is purely imaginary such that the sec- 
ond term in the denominator is always positive and the trans- 
fer function has no poles. If <5 2 w 4 fc 2 _2 (l - e ^ 2ik2d ) < 1, 



ie. for small ku, we can neglect the <5-correction in the de- 
nominator and find an o(l) behavior of the transmission func- 
tion. In the opposite limit of large k\\ we can neglect the one 
in the denominator and find the transfer function decaying 
exponentially with k 2 d. The asymptotic exponential decay 
t- ~ k 2 /(5 2 oj 4 ) cxp(— 2fc||rf) can be used to define the cross- 
over momentum fc max d = — log (ui 2 8/k max ) which has an 
explicit solution in terms of the product-log function, 



k max d = -W(-uj 2 5d) 



■ log 8 



(6) 



for small 8. Since 8 is assumed to be of the order of the 
discretization mesh constant this qualitatively represents the 
logarithmic dependence on the discretization observer in the 
TMM study above. 

In conclusion we investigated the transfer function of the 
discretized perfect lens by means of FDTD and TMM simu- 
lations. The TMM has the advantage of computing the trans- 
fer function directly in (fcii , w)-space as well as eliminating 
the problems associated with the explicit time dependence in 
the FDTD simulations. We argue that the peak observed near 
fc|| max in the FDTD transfer function is due to finite time ar- 
tifacts; it does not exist in the TMM simulations. Further we 
found that the finite discretization mesh acts like an imagi- 
nary deviations from the \i = e = —1 of the PL and leads to 
a cross-over in the transfer function from o(l) to exponential 
decay around a maximum parallel momentum fc|| max limiting 
the attainable super-resolution of the PL. We propose a simple 
qualitative model to describe the impact of the discretization 
in terms of effective thin /i-only and e-only surface layers ex- 
posed by the discretized LH slab which have a thickness <5 that 
is of the order of the discretization mesh constant. k\\ . max d is 
found to depend logarithmically on the mesh constant in qual- 
itative agreement with the TMM simulations. Since virtually 
all simulation solve discretized Maxwell equations they are all 
subject to this restriction. 
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